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ABSTRACT 

The utility of a recently developed analytical micromechanics model for the response of 
metal matrix composites under thermal loading is illustrated by comparison with the results gen- 
erated using the finite-element approach. The model is based on the concentric cylinder assem- 
blage consisting of an arbitrary number of elastic or elastoplastic sublayers with isotropic or 
orthotropic, temperature-dependent properties. The elastoplastic boundary-value problem of an 
arbitrarily layered concentric cylinder is solved using the local/global stiffness matrix formula- 
tion (originally developed for elastic layered media) and Mendelson’s iterative technique of suc- 
cessive elastic solutions. These features of the model facilitate efficient investigation of the 
effects of various microstructural details, such as functionally graded architectures of interfacial 
layers, on the evolution of residual stresses during cool down. The available closed-form expres- 
sions for the field variables can readily be incorporated into an optimization algoritlmi in order 
to efficiently identify optimal configurations of graded interfaces for given applications. Com- 
parison of residual stress distributions after cool down generated using finite-element analysis 
and the present micromechanics model for four composite systems with substantially different 
temperature-dependent elastic, plastic and thermal properties illustrates the efficacy of the 
developed analytical scheme. 


INTRODUCTION 

The concentric cylinder geometry has been employed quite frequently during the past thirty 
years for modeling the response of unidirectional composites. The early models are based on a 
single fiber embedded in a cylindrical shell (Hashin, 1962; Hill, 1965; Mulhem and Rogers, 
1967; Dvorak and Rao, 1976a,b). These models have been used to predict elastic moduli or their 
bounds for general, three-dimensional loading, as well as thermoelastic and inelastic response 
under axisymmetric loading. More recently, generalizations of these models have been proposed 
in order to include the interaction of the composite cylinder with the surrounding medium 
(Dvorak and Bahei-El-Din, 1979; Christensen and Lo, 1979), or to take into account the effect of 
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microstructural details on the overall composite response (Hashin, 1990; Benveniste et al., 1991; 
Sutcu, 1992; Warwick and Clyne, 1991). The microstructural details include the presence of car- 
bon coating around certain types of fibers, the layered morphology of such ceramic fibers as the 
SCS6 SiC fiber used in titanium matrix composites, as well as the presence of an interface or 
interphase layer between the fiber and the matrix. The interfacial region can arise naturally due 
to chemical reaction between the fiber and the matrix (Wawner and Gundel, 1991), or can be 
introduced deliberately in order to minimize residual fabrication stresses (Arnold et al., 1990). 
With the exception of the model proposed by Sutcu which deals with elastic response of uni- 
directional metal matrix composites with layered SiC fibers, the aforementioned models focus 
on specific geometries involving only a few concentric cylinders. 

Recently, an analytical elastoplastic solution to the micromechanics concentric cylinder 
model under axisymmetric thermal loading that allows consideration of arbitrarily layered con- 
figurations has been constructed by Pindera et al. (1992). The solution is based on the 
local/global stiffness matrix formulation originally developed by Bufler (1971) for stress 
analysis of isotropic, elastic layered media, and Mendelson’s method of successive elastic solu- 
tions for elastoplastic problems with Prandtl-Reuss constituents (Mendelson, 1983). Applica- 
tions of the local/global stiffness matrix technique to the elastic stress analysis problems of com- 
posite materials and structures have been presented by Pindera (1991). Derstine and Pindera 
(1988) used the method to develop an approximate solution to the problem of an arbitrarily lam- 
inated composite tube under axisymmetric loading with inelastic graphite/epoxy plies that were 
modeled using the endochronic theory. 

In the local/global stiffness matrix formulation, the solution of a given boundary-value 
problem is reformulated in terms of the interfacial displacements as the basic unknowns in place 
of the eigenvectors associated with the eigenfunctions that satisfy the governing differential 
equations for the given problem. This reformulation entails construction of a local stiffness 
matrix that relates the tractions at the outer boundaries of a layer to the corresponding interfacial 
displacements. By assembling the local stiffness matrices along the main diagonal of the global 
stiffness matrix in an overlapping fashion, the interfacial continuity conditions and the external 
boundary conditions are identically satisfied. In elastoplastic problems, the local stiffness matrix 
equation that describes the response of a layer, and thus the global stiffness matrix for the entire 
assemblage, involves integrals of plastic strain distributions in the individual layers. These 
integrals appear as vectors on the right hand side of the system of equations that ensures satisfac- 
tion of boundary conditions and continuity of interfacial tractions and displacements. Since these 
integrals depend implicitly on the unknown interfacial displacements, the system of equations 
must be solved in an iterative fashion. The solution developed by Pindera et al. (1992) utilizes 
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the classical incremental plasticity model for the inelastic response of the matrix phase in an uni- 
directional metal matrix composite, and employs Mendelson’s method of successive elastic solu- 
tions to solve the system of nonlinear equations. 

The reformulation of the boundary-value problem in terms of the interfacial displacements 
as the basic unknowns using the local/global stiffness matrix approach offers a number of advan- 
tages over the standard formulation. First of all, the construction of the global stiffness matrix 
can be automated very easily, facilitating addition of extra layers without any difficulty. This 
construction also eliminates certain redundant equations that arise from the application of inter- 
facial continuity and boundary condition equations in the standard formulation, resulting in 
nearly 50% reduction in the size of the system of simultaneous equations for a large number of 
layers. Furthermore, as the elements of the stiffness matrices for different types of elastic (isotro- 
pic, transversely isotropic, orthotropic) and inelastic (isotropic) layers have been provided in 
closed form, a given boundary-value problem does not have to be resolved each time a particular 
concentric cylinder assemblage is considered. Different configurations are efficiently handled by 
assembling the global stiffness matrix in an appropriate fashion using the provided local stiffness 
matrices. 

The elastoplastic solution of the multiple concentric cylinder (MCC) under axisymmetric 
thermal loading based on the local/global stiffness matrix formulation is very well suited for 
computer implementation, facilitating efficient parametric studies for developing new composite 
materials with elastoplastic constituents. Since the solution methodology is analytical rather than 
numerical (eg. finite-element or finite-difference schemes), there is no need to generate meshes 
or grids every time the geometrical details of the concentric cylinder are changed. Different con- 
figurations are conveniently handled by changing a few lines in the input file of a computer pro- 
gram that can be executed on a personal computer/work station. It is this feature that makes 
computer implementation of the solution very attractive for use by the materials scientist, 
designer and analyst alike. The construction of the input data file is simple and can be further 
facilitated by a user-friendly interface that allows the user to specify the input parameters in an 
interactive fashion through menu-driven data entry. 

The ability to easily vary the number and arrangement of the concentric cylinders makes 
this method ideal for a wide variety of applications involving axisymmetric loading of advanced 
unidirectional composites. These applications include the reduction of residual stresses in metal 
matrix composites using single or multiple compliant/compensating interfacial layers as sug- 
gested by Arnold and co-workers (1990, 1992) (see also Jansson and Leckie (1992)), and design 
of engineered interfaces for improved performance (Arnold and Wilt, 1992). Further, since the 
solution is analytic, with expressions for the various field quantities provided in closed form, it 
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can readily be incorporated into an optimization algorithm in order to efficiently identify optimal 
interfacial layer configurations or morphologies for given applications. The method also facili- 
tates efficient discretization of functionally graded composites into an arbitrarily large number of 
layers in order to model spatial property variations as accurately as is desired. Thus either 
discretely-graded or continuously-graded interfacial regions can be investigated. 

Thus far, the MCC model has been employed to investigate the effect of: 1.) morphology of 
SCS6 SiC fibers; 2.) architecture of interfacial layers; and 3.) microstructure of the a 2 + P 
ti tani um aluminide matrix, on the evolution of residual stresses in SiC/Ti 3 A1 unidirectional com- 
posites (Pindera et al., 1992; Pindera and Freed, 1992). A recent investigation addresses the 
effectiveness of graded interfacial layers in reducing residual stresses in titanium matrix compo- 
sites subjected to both monotonic and cyclic thermal loading (Williams et al., 1993). An investi- 
gation addressing the problem of optimization of residual stresses in metal matrix composites 
using multiple layers at the fiber/matrix interface is currently in progress. 

A limited comparison between the predictions of the MCC model and the finite-element 
results generated using a commercially available code was presented by Pindera et al. (1992) for 
the problem of a composite cylinder with a single interfacial layer subjected to a temperature 
drop. The purpose of this paper is to present an extensive comparison between the predictions of 
the MCC model and the finite-element results obtained with the commercially-available code 
ABAQUS (1989) for several concentric cylinder configurations with single and multiple interfa- 
cial layers, and different constituents, subjected to thermal loading. The ultimate objective is to 
demonstrate the utility of the developed analytical solution methodology for efficiently and 
accurately investigating the thermoplastic response of MMC’s under axisymmetric thermal load- 
ing in the presence of various microstructural details. 

ANALYTICAL MODEL 

The analytical model is based on a long, cylindrical assemblage of an arbitrary number of 
concentric cylinders or shells perfectly bonded to each other. Figure 1. Each of the cylindrical 
shells is either elastic or inelastic. The elastic shells may be isotropic, transversely isotropic, or 
orthotropic (radially or circumferentially), while the inelastic shells are taken as initially isotro- 
pic and are modelled using the time-independent incremental plasticity with isotropic hardening. 
All the material parameters governing the response of the elastic and inelastic layers arc func- 
tions of temperature. The geometrical model is thus general enough to allow one to model the 
actual morphologies of fiber, interfacial layer and matrix phases in sufficient detail. 

A distribution of displacements and stresses in the individual phases of the concentric com- 
posite cylinder model is sought under the conditions of a spatially uniform temperature change 
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that varies with time. A solution to the outlined elastoplastic boundary-value problem is obtained 
using the displacement formulation. In solving the boundary-value problem, the following nota- 
tion is adopted. The inner solid core is denoted by a subscript or superscript 1 and the outermost 
cylindrical shell by n. The inner radius of the kth shell is denoted by r k _! and the outer radius by 
r k . The traction and displacement components at the inner and outer radii of the kth shell are 
assigned superscripts and respectively. 

For the prescribed axisymmetric loading, the longitudinal, tangential and radial displace- 
ment components u, v and w, referred to the cylindrical coordinate system x-r-8 centered at the 
origin of the concentric cylinder assemblage have the form, 


u = u(x) = Eq x 


w = w(r) 


( 1 ) 


v = 0 


where Eo is the same uniform longitudinal strain for all layers. These displacement components 
yield the following strain components in the cylindrical coordinate system. 


_ du 

£xx_ dx ^ 


e ee = 


w(r) 


£n = 


dw(r) 

dr 


( 2 ) 


with the shear strain components identically zero. Since the strain components are either con- 
stant or functions of only the radial coordinate r, the stress components are at most functions of 
r, and therefore the stress equilibrium equations in cylindrical coordinates reduce to the single 
equation, 

dGrr + Grr ~ Gee = o (3) 

dr r 


The governing differential equation for the radial displacement w(r) in each shell is 
obtained by expressing the stress components On- and Gee in Equation (3) in terms of w(r) and its 
gradient using stress-strain equations and strain-displacement relations given by Equation (2). 
For problems in cylindrical coordinates, the stress-strain equations for an orthotropic material in 
the presence of thermal loading and inelastic effects, and in the absence of shear strains, are 
given by. 
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( 4 ) 




Qa Qi0 Qir 


E xx ®xx(T Tq) 

°ee 

► = 

Cx0 Qe Or 

4 

e ee -eee-aee(T-To) ► 



c„ Cer Crr 


^-ejj-o^CT-To) 


In the above, e xx , e e0 , £n are total strains, e£, e&, e“ arc inelastic strains, and a^T-To), 
ocee(T - T 0 ), ct^CT - T 0 ) are thermal strains, with T 0 denoting a reference temperature and T 
denoting the current temperature. The corresponding constitutive equations for an isotropic or 
transversely isotropic material can be obtained by equating the appropriate stiffness matrix ele- 
ments. 

Given the functional form of the displacement field, the strain-displacement equations and 
the constitutive equations in a cylindrical shell, the governing differential equation for the radial 
displacement w(r) in each shell is obtained from the surviving equilibrium equation in the form. 


Transversely isotropic, elastic layers ( C T a = C^, Qm - C^-, = (y ) 


d 2 w 1 dw _ _w 
dr 2 + r dr r 2 


(5a) 


Qrthotropic, elastic layers 


d 2 w 1 dw 
dr 2 + r dr 


1 Q)0 1 r (Cex Qu-) 

7 c7 w= 7 [ oi 


£0+ X 

i=x,0,r 


(Qi-Cei) 

Or 


0Cii(T-T 0 )] 


(5b) 


Isotropic, inelastic layers 


d 2 w 1 dw 
dr 2 + r dr 


1 1 v 

7 7 2 


(Qi-C ei) 

Cn 


eSi « + 


dr 




i=x,0, 


j ^n- 


(5c) 


where the distribution of the inelastic strains, ejf (r), is assumed to be known at the beginning of 
each thermal load increment. 
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The solution of the above equations is obtained subject to the external boundary condition, 
On-(r n ) = 0, the interfacial displacement and traction continuity conditions, 
w k-i( r k-i) = w k(r k -i), ojT 1 (r k _i) = a^(r k -i), and the longitudinal equiUbrium condition, 
J o xx dA c = 0, where Ac is the cross-sectional area of the concentric cylinder assemblage. 

4 

Using standard techniques, solutions to the governing differential equations are obtained in 
the form. 

Transversely isotropic, elastic layers 


A 2 

w(r) = Air + — 


(6a) 


Orthotropic, elastic layers 


w(r) = Aii^ + A 2 r ^ + 


(Qk Qnt) 

- Cee) 


reo + 


X 

i=x,0,r 


(Cri-Cfli) 
(Or -C ee ) 


Oiir(T-T 0 ) 


(6b) 


Isotropic, inelastic layers 


w(r) = ^Ii 


1 f „ (Cri + Cei) . r f ^ (Q + 


r k _j i=x, 0,r 


Q, 


-ejf (r')r' dx' + i| £ 


r xfc_ii=x,6,r 


Or 


(6c) 


2 

A 2 1 _ Cn ^ / r k— 1 .v 

Ajr+ V — X 7^ e ii ( r k-l M _2 ^ 

r 1 i=x,0,r '-nr ‘ 


where X, = (^r~) 112 and r k _ } < r < r k . 

The above solutions contain unknown coefficients Ai and A* for each layer, as well as the 
unknown, uniform axial strain £q. For the solid core, the constant A 2 vanishes since the radial 
displacement at the center has to vanish. These unknown coefficients are determined from the 
boundary condition, interfacial traction and displacement continuity conditions, and the longitu- 
dinal force equilibrium condition. Application of these conditions yields a system of equations in 
the unknown Ai and A 2 coefficients and the uniform longitudinal strain £q, that is solved 
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iteratively when the inelastic strains are present. An iterative procedure is required because the 
inelastic strains depend implicitly on the coefficients Ai and Ai- 

In order to automate the construction of this system of equations so that any arbitrarily lay- 
ered configuration can easily be considered, we reformulate the problem in terms of the interfa- 
cial radial displacements as the basic unknowns in place of the coefficients A k and A k by using 
the concept of a local stiffness matrix. The local stiffness matrix relates the interfacial tractions 
at the inner and outer radii of the kth layer to the corresponding interfacial radial displacements, 
and is obtained from the solutions to Equations (5a-c) (i.e. Equations (6a-c» in conjunction with 
the constitutive equations and strain-displacement equations. To construct the local stiffness 
matrix for the kth layer, we first express the coefficients A k and A k in terms of the interfacial 
displacements w k( r k-i ) and w^r^) by evaluating the solution for the radial displacement com- 
ponent w(r) at the appropriate locations. These expressions are then used in the equation for the 
radial stress component in the kth layer given in terms of the determined radial displacement 
field. The final step entails an evaluation of the radial stress in the kth layer at the inner and 
outer radii in order to generate the radial tractions at those locations. 

The form of the local stiffness matrix equation for the kth layer in the state of generalized 
plane strain and in the presence of thermal and inelastic effects is 


^rr 

k 


ki2 

k 

w 


* 13 ' 

1" 1 

ff 


[g.V 

K 

1 = 

*21 

1 

CM 

M 

1 

w + 

. .1 

- + < 

'k • 

M 

" eo 

h 

t J 

■ (T-Tq) + < 

W 


The elements k k i , ... , k k 3 of the local stiffness matrix are functions of the geometry and elastic 
material properties of the kth layer (which may vary with temperature). The thermal effects are 
represented by f k and f k , which are functions of the thermal expansion coefficients for the kth 
layer. The plastic effects are represented by g k and g k > which are given in terms of the integrals 
of the plastic strain distribution in the given layer. The elements of the local stiffness matrix and 
the elements of the force vectors appearing in Equation (7) have been provided by Pindera et al. 
(1992) for transversely isotropic and orthotropic, elastic layers, and isotropic, inelastic layers. 

Imposition of continuity of displacements and tractions along the common interfaces, 
together with the boundary condition on the radial stress at r„, and the longitudinal equlibrium 
condition, gives rise to a system of equations in the unknown inteifacial displacements expressed 
in the matrix form below. The first n equations result from the imposition of interfacial con- 
tinuity and boundary conditions, whereas the nth+1 equation ensures that the longitudinal equili- 
brium is satisfied at any cross-section along the cylinder’s axis. 
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k^+kn k ?2 0 • k 23 +k ?3 

k2i k22 +k ll • 


Wi 

w 2 


fl+fl 


g2+gl 

0 k 2 i 

• k^ kS 3 


W n 

£o 

► 4 

f n 2 

E^k 

* J 

>(T-T 0 ) - - 

gS 

E n k 


The elements of the last row, <J>n , <t> 2 2 > Vk. and n k have also been provided by Pindera et al. 
(1992). We observe that the global stiffness matrix is constructed by first superposing the local 
stiffness matrices along the main diagonal in an overlapping fashion, and then adding a column 
and a row to account for the thermal effects and the longitudinal equilibrium condition in the 
case of free thermal expansion/contraction. Under the conditions of plane strain, £<3 vanishes and 
the nth+1 row and column are not added into the global stiffness matrix. It is a simple matter to 
construct a computer algorithm for assembling the global stiffness matrix. 

SOLUTION PROCEDURE 

The system of equations given by Equation ( 8 ) is solved iteratively at each temperature step 
for the specified loading after the manner suggested by Mendelson (1983). The iteration is per- 
formed on the plastic force vector that consists of the elements gj, g 2 and E^k- Th e elements 
gi and g 2 are expressed in terms of the integrals of the plastic strain distributions in the given 
layer that have the form (see Pindera et al. (1992)), 


ik 


1.x 

r k _ii=x,0,r 


(Qj + Cgj) 
Or 


e|?(r')r'dr' 


Tk 

1.x 

r k-i 


(C ri -C &1 ) U i dr / 

Cry r ' 


(9a) 


whereas the element Il k has the form. 


*k 

J E C* i eif(r')r , dr 

r k _!i=x,0,r 


(9b) 


Since the elements of the global stiffness matrix at a given temperature are constant, only one 
inversion of the matrix for each sequence of iterations is required. As the elements fj, f 2 ^ 
££ 2 k of the thermal force vector are also constant at a given temperature, most of the 
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computational effort lies in evaluating the integrals in Equation (9) at each iteration. The algo- 
rithm for the iterative procedure is given in the sequel. 

For the given temperature increment, the plastic strain distribution in each layer is 
expressed in terms of the distribution at the preceding temperature plus an increment that results 
from the imposed temperature change. 


efj(r) I current efj(r)l previous + de{j(r) (10) 

The plastic strain increment is derived from the von Mises yield condition which, in the presence 
of temperature-dependent elastoplastic properties of the matrix phase, has the form 

F = yC-jO-j - jo 2 (f,T) = 0 (11) 


where g is the effective yield stress, which is a function of both the effective plastic strain e and 
temperature. The plastic strain increment is thus, 



-^-dX = GijdX 


( 12 ) 


where the proportionality constant dX is obtained from the consistency condition for plastic load- 
ing in the form 


Gij D ijki(d£ki ~ <fc£) + ( G ijCij )^ T 


dX = 


« ,J 3 3T 


4-2 do , / 

— G — p- + GjjDijkiGid 

9 de 


(13) 


where dX > 0 for plastic loading, dX < 0 for neutral loading or unloading, and de* are the ther- 
mal strain increments given by. 


de* = [<Xij(T) + 


aa,j(T) 

dT 


(T-T 0 )]dT 


(14) 
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dD* 

D ijkl are the elastic stiffness elements, and C-j = -^jr-£ki- 10 * e present investigation, the elas- 
toplastic stress-strain response of the matrix is taken to be bilinear so that the slope of the effec- 
tive stress-plastic strain curve, — , is constant at a given temperature. Since the incremental 

d£p 

theory of plasticity is used to calculate the plastic strain increment at each point along the radial 
coordinate, the inelastic analysis is valid for both loading and unloading in the plastic phases. 

The plastic strain distribution in each layer is determined by calculating the plastic strains 
at twenty-one stations after updating the plastic strains at these locations using Equation (10). 
The current values for the plastic strains at these stations are then used in determining the 
integrals given in Equation (9), and thus the elements of the plastic force vector in Equation (8). 
Updated values of the interfacial displacements are then obtained from Equation (8). With a 
knowledge of the interfacial displacements and the axial strain Eq, the coefficients A^andA^in 
each layer can be obtained, producing solutions for the radial displacement (r) from which 
radial and tangential total and plastic strains, and the corresponding stresses, can be obtained. 
These are then used to obtain new approximations for the plastic strain increments. The iterative 
process is terminated when the differences between two successive sets of plastic strain incre- 
ments are less than some prescribed value. The above procedure is described in detail by Men- 
delson (1983). 

NUMERICAL RESULTS AND DISCUSSION 

In order to demonstrate the accuracy, efficiency and efficacy of the outlined solution tech- 
nique, we investigate the evolution of residual stresses in unidirectional metal matrix composites 
with and without interfacial layers for four material systems with different temperature depen- 
dence of the constituents’ behavior. The predictions of the MCC model subsequently are com- 
pared with the results generated using the commercially available finite-element code ABAQUS. 

The considered material systems are state-of-the-art metal matrix composites under con- 
sideration for advanced aerospace applications. The four material systems chosen for the correla- 
tion study are: SiC/Ti 3 Al, Al 2 0 3 /NiAl, Al 2 0 3 /NiCrAlYl, and Al 2 0 3 /NiCrAlY2. The properties 
of the indvidual constituents for these composites are listed in Table 1 (SCS6 fiber and Ti-24A1- 
llNb matrix) and Table 2 (A1 2 0 3 fiber and NiAl, NiCrAlYl, and NiCrAlY2 matrices). For the 
SiC/Ti 3 Al system, residual stresses were calculated for cool down from 1500°F to 75°F. In the 
case of the remaining material systems, residual stresses were calculated for cool down from 
2300°F to 70°F. 
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For the SiC/ri 3 A1 system, three configurations were considered, namely; a fiber embedded 
in a homogeneous matrix, a fiber embedded in a matrix with a one-layer interface, and a fiber 
embedded in a matrix with a three-layer interface. In all cases, the outer radius of the composite 
cylinder was normalized to 1.0, with the normalized fiber radius of 0.6320 producing a fiber 
volume fraction of 0.40. For the configurations involving one- and three-layer interfaces, the 
outer radius of the interfacial region was 0.6952, resulting in an interfacial volume fraction of 
0.08. The mechanical properties of the interface layers (with the exception of the Poisson’s ratio) 
were taken to be one half of the matrix properties at each temperature. The thermal expansion 
coefficient of the single interfacial layer configuration was two times that of the matrix, while in 
the case of the three-layer interface configuration the thermal expansion coefficient was two, 
two and one half, and three times that of the matrix phase starting from the innermost layer. For 
the AI 2 O 3 fiber systems, only one configuration was considered, namely; a fiber embedded in a 
matrix without an interface having the same geometry as the corresponding SiC/T^Al compo- 
site. 

The MCC calculations were performed using temperature increments of - 1.0°F. Conver- 
gence of plastic strain increments at the various radial locations typically did not require more 
than six iterations at each temperature increment. As an additional check, values of the effective 
stress calculated from the effective stress-plastic strain curve at various radial locations were 
compared with values of the effective stress based on the obtained stress components at these 
locations. Typically, differences were a fraction of a percent. 

For the finite-element calculations, three quarter-symmetry meshes were constructed that 
consisted of triangular (CGPE 8 ) and quadrilateral (CGPE10), generalized plane strain elements, 
Figure 2. The coarse mesh (Figure 2a) consisted of 88 elements, the medium mesh (Figure 2b) 
consisted of 296 elements, and the fine mesh (Figure 2c) consisted of 580 elements. Each of the 
three meshes has provision for inclusion of a single interfacial layer with four elements across 
the thickness of the interfacial region. These meshes were used to study the convergence and 
accuracy of the finite-element results and to compare the CPU times with the corresponding run 
times for the MCC calculations. Selected comparison of the axial, circumferential and radial 
stress components and the radial displacement at several radial locations in the fiber, interfacial 
layer and matrix obtained with the present analytical solution and the finite-element solution is 
presented in Table 3 for the SiC/Ti 3 A1 material system with one-layer interface. The results indi- 
cate that even with the coarse mesh the differences between the MMC and finite-element results 
are quite small. The fine mesh provides better approximation of the external zero-traction boun- 
dary condition and in some cases better satisfaction of the interfacial continuity conditions. For 
this reason, it was used to generate the remainder of the results for the SiC/T^Al system with 
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one-layer interface, as well as for the configuration without an interface. The same mesh was 
also used to generate results for the AI 2 O 3 fiber systems. For the SiC/T^Al system with a 
three-layer interface, a new mesh was constructed with 708 elements that was based on the fine, 
580-element mesh. Figure 3. This was done by subdividing the interfacial region into three 
layers while maintaining four elements through the thickness of each of the three layers. Table 4 
presents a comparison of CPU times for the cool down of the SiC/T^Al composite with one- 
layer and three-layer interfaces from 1500°F to 75°F simulated with the four finite-element 
meshes and the MCC model. The calculation were performed on an IBM RS-6000 machine. The 
efficiency of the MCC model is clearly evident. 

It should be noted that the comparison of the CPU times presented above is perhaps 
misleading because the analytical model is one-dimensional, whereas the finite-element calcula- 
tions were performed using a two-dimensional model. It would certainly have been more 
appropriate to generate the finite-element results using axisymmetric, generalized plane strain 
(line) elements. Such elements, however, are presently not available in the existing ABAQUS 
code. Nevertheless, the comparison presented above does give an indication of the model s effi- 
ciency. Further, while the finely-refined meshes used in the present investigation are probably 
not necessary for most applications, and may even be considered an overkill , the authors 
nevertheless are of the opinion that in order to meet the intended objectives of this investigation, 
such refinement was desirable. This is because the outlined analytical solution is, in principle, an 
exact solution that satisfies, in addition to the field equations, all the external boundary condi- 
tions and internal continuity conditions provided that the integration of the plastic strain distribu- 
tions is performed accurately. In the present case, the interfacial region was subdivided into 
three layers of equal thickness and the matrix phase into six, and the plastic strain distributions 
were calculated at twenty-one points within each layer as previously mentioned. The accuracy 
with which the integrals of the plastic strain distributions given in Equation (9) were calculated 
in generating the analytical solution, therefore, was quite high. Consequently, in order to obtain 
a comparable level of resolution with the finite-element approach and thus to render the com- 
parison between the two solutions meaningful, the finely-refined mesh configuration was 
chosen. 

SiC/Ti^Al system 

Figure 4a-c presents the axial, circumferential and radial stress distributions, respectively, 
in the SiC/Ti 3 A1 unidirectional composite without an interface, at four temperatures between the 
stress-free temperature of 1500°F and the final cool-down temperature of 75°F. The analytical 
predictions are indicated by the solid lines whereas the finite-element results are given by 
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different symbols superposed on the solid lines corresponding to the four different temperatures. 
The stress distributions at 875°F are elastic, characterized by piecewise uniform values of the 
axial stress throughout the entire region, and inversely quadratic dependence of the circumferen- 
tial and radial stresses in the matrix phase. These baseline distributions can be compared with 
the distributions at the remaining temperatures which exhibit significant influence of plasticity at 
the fiber/matrix interface as discussed by Arnold et al. (1990). Comparing the distributions at the 
different temperatures allows one to delineate the elastic-plastic boundaries in the matrix phase. 
Clearly, the correlation between the analytical solution to the MCC model and the finite-element 
results is excellent. 

The three stress distributions for the SiC/T^Al composite with a one-layer and three-layer 
interface at 75°F after cool down from 1500°F are compared in Figure 5a-c. Again, the correla- 
tion between the analytical solution to the MCC model and the finite-element results for the two 
configurations is unquestionably excellent. Comparing the stress distributions for the one-layer 
and three-layer interface configurations one observes the effect of grading the thermal expansion 
coefficient in the interfacial region on the stress fields in the matrix phase as well as the interfa- 
cial region itself. The use of compensating layers in reducing fabrication-induced residual 
stresses in unidirectional metal matrix composites in the presence of matrix plasticity has been 
discussed by Arnold et al. (1990, 1992), as well as Pindera et al. (1992), and is currently the sub- 
ject of continuing investigation. 

AliOt based systems 

Figure 6 a-c presents the distributions of the axial, circumferential and radial stresses at 
70°F after cool down from 2300°F for the three AI 2 O 3 fiber-based material systems. The two 
NiCrAlY matrices, designated by numbers 1 and 2, differ only in the temperature dependence of 
the Young’s moduli and yield stresses. The Young’s modulus of the NiCrAlY 1 matrix is 
between 38% and 55% higher than that of NiCrAlY2 in almost the entire temperature range 
(with the exception of 2400°F), while the range in which the yield stress of NiCrAlY 1 exceeds 
that of NiCrAlY2 varies between 25% and 155%. Comparing the thermoelastoplastic properties 
of the NiAl matrix on the one hand, and the two NiCrAlY matrices on the other, one observes: 
initially similar magnitudes in the thermal expansion coefficient but a steeper rise with increas- 
ing temperature for the NiCrAlY matrices; similar magnitudes of the Young’s modulus in the 
lower temperature range for the NiAl and NiCrAlY matrices, but a greater decline with increas- 
ing temperature for the NiCrAlY matrix; significantly lower yield stress in the NiAl matrix up to 
1800°F; and the hardening slope of the NiAl matrix that is five times greater than that of the two 
NiCrAlY matrices. 
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The results in Figure 6a-c indicate that the NiCrAlYl matrix composite produces the 
highest axial, circumferential and radial residual stresses, followed by the NiCrAlY2 matrix 
composite and then finally the NiAl matrix composite. The differences in the three residual 
stress distributions between the two NiCrAlY matrix composites and the NiAl matrix composite 
are significant. The stress distributions in the three material systems exhibit significant plasticity 
effects, with the NiAl matrix system showing the most uniform profiles in the matrix region. 
This is not surprising in view of the low values of yield stress in comparison to the yield stress of 
the two NiCrAlY matrices. Comparing the analytical solution of the MCC model and the finite- 
element results, one again observes excellent correlation. 

CONCLUSIONS 

The extensive comparison between the results of the analytical micromechanics model 
based on the multiple concentric cylinder geometry and the finite-element analysis for several 
different material systems exhibiting different temperature dependence of the constituents’ pro- 
perties and different configurations, clearly demonstrates the accuracy and the power of the 
developed analytical approach for modeling thermoplastic response of unidirectional metal 
matrix composites under axisymmetric loading. The pronounced plasticity effects observed in 
the response of advanced metal matrix composites currently under development for elevated 
temperature applications necessarily require incorporation of inelastic constitutive models for 
the matrix and/or fiber phases into the micromechanical stress analysis. Numerical procedures 
such as the finite-element or finite-difference analyses are approximate, geometry-specific and 
often require large computation times. The current method is an extension of the local/global 
stiffness matrix formulation originally developed for stress analysis of layered, elastic media, to 
axisymmetric elastoplastic problems and utilizes Mendelson’s method of successive elastic solu- 
tions. The major advantage offered by the outlined method is the ease with which different 
geometric configurations involving arbitrarily layered concentric cylinders are handled. Since 
the method is analytic, this is easily accomplished by changing a few lines in the input data file 
of a computer code. This, in turn, facilitates efficient parametric studies in the course of develop- 
ing new composite materials, as well as in investigating the effects of various microstructural 
details on the local and global response. The demonstrated predictive capability and efficiency 
of the analytic MCC model sets the stage for further investigations involving unidirectional 
metal matrix composites with functionally-graded constituents. 
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Table 1. Material properties of SCS6 SiC fiber and titanium matrix (Arnold et al. (1990)). 


Material properties 

75°F 

392°F 

797°F 

1112°F 

1202°F 

1500°F 

SiC fiber 

a(KT 6 /°F) 

1.96 

2.01 

2.15 

2.33 

2.38 

2.50 

E (Msi) 

58.00 

58.00 

58.00 

58.00 

58.00 

58.00 

V 

0.25 

0.25 

0.25 

0.25 

0.25 

0.25 

H-24A1-1 INb matrix 

a (10 -6 / °F) 

5.00 

5.20 

5.70 

5.85 

5.90 

6.15 

E (Msi) 

16.00 

14.50 

11.00 

12.50 

9.89 

6.20 

V 

0.26 

0.26 

0.26 

0.26 

0.26 

0.26 

Gy (ksi) 
H (Msi) 

53.89 

59.00 

53.70 

42.20 

39.10 

24.00 

3.33 

0.44 

0.32 

0.19 

0.097 

0.00 


Note: In the table above and the tables that follow, a is the instantaneous thermal expansion 
coefficient; E is the Young’s modulus, v is the Poisson’s ratio; c y is the yield stress; and H 
is the hardening slope based on a bilinear representation of the elastic-plastic stress-strain 
response. 
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Table 2. Material properties of A1 2 0 3 fiber and NiAl (Bowman, 1992), NiCrAlYl, and NiCrAlY2 

(Barkalow et al., 1985) matrices. 


Properties 

70°F 

800°F 

1200°F 

1400°F 

1600°F 

1800°F 

2000°F 

2400°F 

A1 2 0? fiber 

a (10 -6 / °F) 

3.35 

4.86 

5.42 

5.63 

5.79 

5.91 

5.98 

5.97 

E (Msi) 

61.10 

59.83 

59.60 

59.33 

59.07 

58.80 

53.13 

50.30 

V 

0.20 

0.20 

0.20 

0.20 

0.20 

0.20 

0.20 

0.20 

NiAl matrix 

aaO -6 /^) 

7.12 

8.24 

8.76 

9.01 

9.23 

9.45 

9.69 

9.95 

E (Msi) 

28.00 

25.59 

24.27 

23.61 

22.95 

22.29 

19.29 

10.14 

v 

0.32 

0.32 

0.32 

0.32 

0.32 

0.32 

0.32 

0.32 

o y (ksi) 

45.70 

22.50 

16.06 

13.07 

10.62 

7.79 

5.30 

3.39 

H(Msi) 

0.50 

0.50 

0.50 

0.50 

0.50 

0.50 

0.50 

0.50 

NiCrAlYl* matrix 

a(KT 6 /°F) 

5.44 

7.76 

9.03 

9.67 

10.30 

10.94 

11.57 

12.84 

E (Msi) 

27.00 

25.30 

21.40 

17.70 

12.80 

7.60 

3.10 

0.10 

v 

0.30 

0.30 

0.30 

0.30 

0.30 

0.30 

0.30 

0.30 

o y (ksi) 

150.00 

146.00 

135.80 

101.20 

74.90 

21.90 

5.00 

1.00 

H(Msi) 

0.10 

0.10 

0.10 

0.10 

0.10 

0.10 

0.10 

0.10 

NiCrAlY2** matrix 

a (10 -6 / °F) 

5.44 

7.76 

9.03 

9.67 

10.30 

10.94 

11.57 

12.84 

E (Msi) 

19.50 

16.80 

14.70 

11.60 

8.90 

5.20 

2.10 

0.10 

v 

0.30 

0.30 

0.30 

0.30 

0.30 

0.30 

0.30 

0.30 

Oy (ksi) 

117.80 

104.10 

99.40 

71.00 

29.40 

11.00 

3.80 

0.70 

H(Msi) 

0.10 

0.10 

0.10 

0.10 

0.10 

0.10 

0.10 

0.10 


* NiCrAlYl matrix chemical composition: Ni-10.2Cr-9.3Al-6.0Ta-022Hf-0.43Y 

** NiCrAlY2 matrix chemical composition: Ni-17.8Cr-12YAT0.61Y-0.73Hf-20.8Co 
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Table 3. Comparison of stresses and radial displacements generated by the present model 

and ABAQUS. 


Radial location 

(psi) 

<5»e 

(psi) 

On 

(psi) 

w(r) 

(xlO 2 in/in) 

D.6320 (fiber) 





ABAQUS (88 el.) 

-83.855 

-19,579 

-19,580 

-0.19012 

ABAQUS (2% el.) 

-83,856 

-19,580 

-19,580 

-0.19012 

ABAQUS (580 el.) 

-83,854 

-19,579 

-19,577 

-0.19012 

Present model 

-83,840 

-19,590 

-19,590 

-0.19010 

Q.6320 (r inn „ interfacial layer) 




ABAQUS (88 el.) 

38,960 

41,897 

-19,557 

-0.19012 

ABAQUS (296 el.) 

38,967 

41,902 

-19,547 

-0.19012 

ABAQUS (580 el.) 

38,973 

41,908 

-19,539 

-0.19012 

Present model 

38,900 

41,870 

-19,590 

-0.19010 

0.6952 (r out , r interfacial layer) 




ABAQUS (88 el.) 

41,792 

39,235 

-14,085 

-0.41143 

ABAQUS (296 el.) 

41,743 

39,190 

-14,141 

-0.41145 

ABAQUS (580 el.) 

41,813 

39,255 

-14,056 

-0.41141 

Present model 

41,720 

39,200 

-14,130 

-0.41150 

0.6952 (matrix) 





ABAQUS (88 el.) 

50,627 

36,758 

-13,831 

-0.41143 

ABAQUS (296 el.) 

50,507 

36,447 

-14,053 

-0.41145 

ABAQUS (580 el.) 

50,445 

36,541 

-14,087 

-0.41141 

Present model 

50,320 

36,640 

-14,130 

-0.41150 

1.0 (matrix) 





ABAQUS (88 el.) 

62,783 

28,045 

97 

-0.72439 

ABAQUS (296 el.) 

62,778 

28,038 

94 

-0.72439 

ABAQUS (580 el.) 

62,714 

28,004 

14 

-0.72439 

Present model 

62,700 

27,990 

0 

-0.72450 
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Table 4. Comparison of CPU times between the present model and ABAQUS for cooling down 
of a SiC/Ti 3 A1 unidirectional composite from 1500°F to 75°F. 


CPU time (seconds) 

1 interfacial layer 


Present model (189 integration stations) 

12.76 

ABAQUS (88-element mesh) 

87.24 

ABAQUS (296-element mesh) 

237.57 

ABAQUS (580-element mesh) 

511.00 

3 interfacial layers 


Present model (189 integration stations) 

12.76 

ABAQUS (708-element mesh) 

714.41 


X 








B.) 296 Element Mesh 



Figure 2 - Finite-element mesh configurations for a concentric cylinder with one interfacial 
layer: a.) 88 element mesh; b.) 296 element mesh; c.) 580 element mesh. 
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C.) 580 Element Mesh 
Figure 2 - Concluded. 



Figure 3 - Finite-element mesh configuration with 708 elements for a concentric cylinder with 
three interfacial layers. 
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matrix 



Figure 4 - Comparison between finite-element and present solution stress distributions in a 
SiC/Ti-24Al-llNb composite without an interface at T = 875°F, 475°F, 275°F and 75°F after 
cool down from 1500°F: a.) axial stress; b.) hoop stress; c.) radial stress. 











Figure 5 - Comparison between finite-element and present solution stress distributions in a 
SiC/Ti-24Al-llNb composite with an one-layer and three-layer interface at T = 75°F after cool 
down from 1500°F: a.) axial stress; b.) hoop stress; c.) radial stress. 
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fiber 


matrix 



A.) 


Figure 6 - Comparison between finite-element and present solution stress distributions 
AI2O3/N1AI, Al203/NiCrAlYl and Al2C>3/NiCrAlY2 composite without an interface at 
70 °F after cool down from 2300 °F: a.) axial stress; b.) hoop stress; c.) radial stress. 
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Figure 6 - Continued. 
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